A Note on the Importance of Weak Convergence Rates for SPDE Approximations in Multilevel Monte Carlo Schemes

It is a well-known rule of thumb that approximations of stochastic partial differential equations have essentially twice the order of weak convergence compared to the corresponding order of strong convergence. This is already known for many approximations of stochastic (ordinary) differential equations while it is recent research for stochastic partial differential equations. In this note it is shown how the availability of weak convergence results influences the number of samples in multilevel Monte Carlo schemes and therefore reduces the computational complexity of these schemes for a given accuracy of the approximations.Comment: 16 pages, 3 figures, updated to version published in the Proceedings of MCQMC1

Regression-based variance reduction approach for strong approximation schemes

In this paper we present a novel approach towards variance reduction for discretised diffusion processes. The proposed approach involves specially constructed control variates and allows for a significant reduction in the variance for the terminal functionals. In this way the complexity order of the standard Monte Carlo algorithm ($\varepsilon^{-3}$) can be reduced down to $\varepsilon^{-2}\sqrt{\left|\log(\varepsilon)\right|}$ in case of the Euler scheme with $\varepsilon$ being the precision to be achieved. These theoretical results are illustrated by several numerical examples.Comment: arXiv admin note: text overlap with arXiv:1510.0314

Discrete adjoint approximations with shocks

This paper is concerned with the formulation and discretisation of adjoint equations when there are shocks in the underlying solution to the original nonlinear hyperbolic p.d.e. For the model problem of a scalar unsteady one-dimensional p.d.e. with a convex flux function, it is shown that the analytic formulation of the adjoint equations requires the imposition of an interior boundary condition along any shock. A 'discrete adjoint' discretisation is defined by requiring the adjoint equations to give the same value for the linearised functional as a linearisation of the original nonlinear discretisation. It is demonstrated that convergence requires increasing numerical smoothing of any shocks. Without this, any consistent discretisation of the adjoint equations without the inclusion of the shock boundary condition may yield incorrect values for the adjoint solution

Blocked All-Pairs Shortest Paths Algorithm on Intel Xeon Phi KNL Processor: A Case Study

Manycores are consolidating in HPC community as a way of improving performance while keeping power efficiency. Knights Landing is the recently released second generation of Intel Xeon Phi architecture. While optimizing applications on CPUs, GPUs and first Xeon Phi's has been largely studied in the last years, the new features in Knights Landing processors require the revision of programming and optimization techniques for these devices. In this work, we selected the Floyd-Warshall algorithm as a representative case study of graph and memory-bound applications. Starting from the default serial version, we show how data, thread and compiler level optimizations help the parallel implementation to reach 338 GFLOPS.Comment: Computer Science - CACIC 2017. Springer Communications in Computer and Information Science, vol 79

A Dimension-Adaptive Multi-Index Monte Carlo Method Applied to a Model of a Heat Exchanger

We present an adaptive version of the Multi-Index Monte Carlo method, introduced by Haji-Ali, Nobile and Tempone (2016), for simulating PDEs with coefficients that are random fields. A classical technique for sampling from these random fields is the Karhunen-Lo\`eve expansion. Our adaptive algorithm is based on the adaptive algorithm used in sparse grid cubature as introduced by Gerstner and Griebel (2003), and automatically chooses the number of terms needed in this expansion, as well as the required spatial discretizations of the PDE model. We apply the method to a simplified model of a heat exchanger with random insulator material, where the stochastic characteristics are modeled as a lognormal random field, and we show consistent computational savings

Eigen-AD: Algorithmic Differentiation of the Eigen Library

In this work we present useful techniques and possible enhancements when applying an Algorithmic Differentiation (AD) tool to the linear algebra library Eigen using our in-house AD by overloading (AD-O) tool dco/c++ as a case study. After outlining performance and feasibility issues when calculating derivatives for the official Eigen release, we propose Eigen-AD, which enables different optimization options for an AD-O tool by providing add-on modules for Eigen. The range of features includes a better handling of expression templates for general performance improvements, as well as implementations of symbolically derived expressions for calculating derivatives of certain core operations. The software design allows an AD-O tool to provide specializations to automatically include symbolic operations and thereby keep the look and feel of plain AD by overloading. As a showcase, dco/c++ is provided with such a module and its significant performance improvements are validated by benchmarks.Comment: Updated with accepted version for ICCS 2020 conference proceedings. The final authenticated publication is available online at https://doi.org/10.1007/978-3-030-50371-0_51. See v1 for the original, extended preprint. 14 pages, 7 figure

Construction of a Mean Square Error Adaptive Euler--Maruyama Method with Applications in Multilevel Monte Carlo

A formal mean square error expansion (MSE) is derived for Euler--Maruyama numerical solutions of stochastic differential equations (SDE). The error expansion is used to construct a pathwise a posteriori adaptive time stepping Euler--Maruyama method for numerical solutions of SDE, and the resulting method is incorporated into a multilevel Monte Carlo (MLMC) method for weak approximations of SDE. This gives an efficient MSE adaptive MLMC method for handling a number of low-regularity approximation problems. In low-regularity numerical example problems, the developed adaptive MLMC method is shown to outperform the uniform time stepping MLMC method by orders of magnitude, producing output whose error with high probability is bounded by TOL>0 at the near-optimal MLMC cost rate O(TOL^{-2}log(TOL)^4).Comment: 43 pages, 12 figure

Rapid covariance-based sampling of linear SPDE approximations in the multilevel Monte Carlo method

The efficient simulation of the mean value of a non-linear functional of the solution to a linear stochastic partial differential equation (SPDE) with additive Gaussian noise is considered. A Galerkin finite element method is employed along with an implicit Euler scheme to arrive at a fully discrete approximation of the mild solution to the equation. A scheme is presented to compute the covariance of this approximation, which allows for rapid sampling in a Monte Carlo method. This is then extended to a multilevel Monte Carlo method, for which a scheme to compute the cross-covariance between the approximations at different levels is presented. In contrast to traditional path-based methods it is not assumed that the Galerkin subspaces at these levels are nested. The computational complexities of the presented schemes are compared to traditional methods and simulations confirm that, under suitable assumptions, the costs of the new schemes are significantly lower.Comment: 18 pages, 5 figures; numerical simulations revised, implementation section added; To appear in Monte Carlo and Quasi-Monte Carlo Methods - MCQMC, Rennes, France, July 201

The VOLNA-OP2 tsunami code (version 1.5)

In this paper, we present the VOLNA-OP2 tsunami model and implementation; a finite-volume nonlinear shallow-water equation (NSWE) solver built on the OP2 domain-specific language (DSL) for unstructured mesh computations. VOLNA-OP2 is unique among tsunami solvers in its support for several high-performance computing platforms: central processing units (CPUs), the Intel Xeon Phi, and graphics processing units (GPUs). This is achieved in a way that the scientific code is kept separate from various parallel implementations, enabling easy maintainability. It has already been used in production for several years; here we discuss how it can be integrated into various workflows, such as a statistical emulator. The scalability of the code is demonstrated on three supercomputers, built with classical Xeon CPUs, the Intel Xeon Phi, and NVIDIA P100 GPUs. VOLNA-OP2 shows an ability to deliver productivity as well as performance and portability to its users across a number of platforms